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ABSTRACT 


An approximate method is presented for analyzing restrained 

shrinkage stresses in ablation heat shields during space flight 

cooling. An alternating direction iterative method is used 'to obtain 

numerical solutions of a two-dimensional mixed boundary- value problem 

for rectangular elastic plates. Plane strain theory is assumed with . 

the heat-.shield length and thicknesses as controlling, dimensions and 

* 

the finite difference method’ is used to integrate the .governing 5 • 

differential equations for displacements and -stresses.. .The calculated 

stresses are shown to agree with exact results well ! away’ from -the heat- 

V • *2 * i 

shield edges. At the restrained-free edge corner of the heat shield, 

1 * 

the calculated stresses do not agree with the prescribed boundary 
stresses. It is concluded that calculated stresses are valid for 
regions other than the corner vicinity, but are not valid near the 
corner. 

Experimental results from cold soak tests of ablation heat- 
shield models are presented. These tests indicate that the location 
of maximum restrained shrinkage stress depends on the ablation material 
thickness. Calculated stresses are qualitatively compared with the 
test data.. The calculations and the test data both indicate that heat 
shields constructed with unreinforced ablation materials are susceptible 
to edge cracking when the ablation material thickness exceeds about 
0. 5 inch. 
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CHAPTER X 


INTRODUCTION 
1.1 Research Problem 

The effects of restrained shrinkage due to space flight cooling 
must be investigated prior to the. materials selection and design of 
ablation heat shields for proposed multimission .manned spacecraft. The 
multimission spacecraft concept involves a single vehicle structure 
which can be adapted for several space missions. The adaptation 
required for these spacecraft is a refurbishable heat- shielding system 
which will consist of ablation material panels with provision for 
removing and replacing the entire system. This thesis presents an 
approximate thermal- stress analysis method which is applicable to 
refurbishable heat shields. The method of analysis, programed for 
solution with a digital computer, provides a means to investigate how 
restrained shrinkage stresses in ablation materials are affected by 
temperature distribution, ablator strength properties, and ablator 
thickness. 

Studies of refurbishable heat shields (Ref. 1) usually consider 
rectangular planform panels of ablation material bonded to substructure 
panels constructed of conventional materials such as steel or fiber- 
glass. Figure 1 shows the arrangement of such panels on a conceptual 
spacecraft and Figure 2 shows some details of an individual panel. 

The ablation materials considered for these panels are specialized 
low-density-plastic composites which have low ultimate tensile 








Figure 2.— Model of ablation, heat-shield panel 
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strengths and high expansion coefficients. Since the heat- shield 
panels on the multimission vehicle will he continuously exposed to the 
space environment, the panels can experience severe cooling due to 
radiation heat transfer. This cooling will produce tension stresses 
in the ablation material because of differential thermal contraction 
between the ablation material and the substrate. The combination of 
tension stress and low ultimate strength can cause the ablation 
material to crack. From a design standpoint, cracks in the ablation 
material prior to reentry are unacceptable because the cracks com- 
promise the structural integrity of the heat shield. Therefore, a 
technique for predicting' restrained shrinkage stresses must be used 
during heat-shield design. Experiments conducted during this thesis 
investigation indicate that the location of maximum shrinkage stress 
varies with thickness and can be greatest at the edges of the heat 
shield. These data indicate the inadequacy of simplified shrinkage 
analyses which neglect edge effects. 

In this thesis, the heat shield is approximated as a flat 
rectangular slab of ablation material bonded oh one side to a stiff 
substructure. Plane strain theory is used with the heat- shield length 
and thickness as the controlling dimensions. The finite difference 
method is used to integrate the governing differential equations for 
displacements and stresses. Temperature gradients through the ablation 
material thickness and variations of elastic constants with temperature 
are permitted in the analysis. Thermal- stress calculations are com- 
pared with exact solutions and are qualitatively compared with 
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experimental data obtained during this investigation. The restrained 
shrinkage analysis involves solving a mixed boundary problem in that 
the bonded surface of the ablation material has prescribed displace- 
ments >11116 all other surfaces have prescribed boundary stresses. 
Experimental and analytical studies of restrained shrinkage or 
expansion have been made by other persons. Some of these studies are 
discussed in the following section. With but one exception, the 
analytical methods used in these studies differ from the thesis method. 

1.2 Historical Background of Research Problem 
Thermal stresses in a rectangular plate, clamped along one edge, 
are analyzed in Reference 2. In this reference, plane stress theory is 
used and the thermal- stress problem with a uniform temperature change 
is converted to a problem of specified fictitious edge stresses. The 
integral representing the strain energy is expressed in. terms, of normal 
and shear stresses and minimized using the calculus of variations. The 
analytical results show a normal stress concentration factor of 10 at. 
the corner junction of the free and clamped edges. The specified 

' i 

interface displacements are, however, not satisfied in this analysis. 

. The analyses of a rectangular plate with a clamped edge reported 
in References 5 and 4 assume plane stress theory and uniform tempera- 
ture. In these analyses, the shear and normal stresses at the 'damped 
edge are expressed in terms of Fourier series satisfying the differ- 
ential equations of equilibrium. Numerical results given in both 
references show infinite normal and shear stresses at the clamped-free 
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edge corner. Reference 4 also presents calculated stresses for plate 
regions other than the clamped edge. The author states that these 
calculated stresses are valid since the corner stress singularity is 
local and confined to the corner vicinity. 

Photoelastic studies of restrained shrinkage are reported in 
References 5? 6, and r J. The problems associated with bonding solid 
rocket propellant grains to rocket motor cases are the bases for these 
studies. The experimental data indicate stress concentrations at the 
clamped- free edge corners. In Reference 5, the author states that 
photoelastic measurements up to 0.001 inch from the corner showed 
rapidly increasing normal tension stresses with approach to the corner. 
The author discusses, but does not reconcile the discrepancy between 
his data and the requirements for compression normal stresses to 
satisfy the boundary conditions for plane strain. 

The analysis of Reference 8 is a plane stress analysis of a 
rectangular plate clamped along one edge and subjected to a uniform 
temperature change. The mixed boundary-value problem is solved for the 
two-dimensional displacement equations of equilibrium using a finite 
difference method. The numerical results show large, but finite values 
for normal and shear stresses at the clamped-free edge corner. The 
author states that the results for other plate locations agree 
qualitatively with the results of Reference ^ and that the precise 
behavior of the solution at the corner is not known. 



CHAPTER II 


ANALYTICAL MODEL 

An analytical model was established to calculate the restrained 
shrinkage stresses in a flat rectangular slab of ablation material 
bonded to a substructure. Figure 3 shows the model configuration and 
identifies the rectangular Cartesian coordinate directions. The 
substructure of an actual heat shield would consist of structural 
elements such as honeycomb sandwich (see Fig. 2) or face sheets and 
stringers, but for this analysis the substructure is simply a 
restraining foundation for the ablation material. 

2.1 As sumptions 

The assumptions used in the analysis are discussed in the 
following paragraphs. 

As shown in Figure 3, the x-y plane is the bonded interface 
between the ablation material and the substructure. The analysis is 
applied to the x-z plane, indicated by the heavy dashed line in Figure 3; 
whose bounded edge coincides with the line -b < x < + b, z = 0. The 
two-dimensional displacement vector of the plane has the components 
u(x,z) in the x-direction and w(x,z) in the z-direction. The two- 
dimensional approximation simplifies the analysis, and it is reasonable 
from consideration of the heat-shield dimensions. The typical x and 
y dimensions for the ablation material on an actual heat shield would 
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be between one and two orders of magnitude greater than the thickness. 
Therefore, if the heat-shield comer regions are excluded, a large 
area remains for which a two-dimensional approximation is applicable. 

It is assumed that the ablation material is fixed to the 
substructure at z = 0 with a bonding material having negligible 
thickness and that the composite panel is free to expand or contract, 
but restrained against curvature. Displacements at the bonded 
interface (z = 0) are determined from the substructure temperature 
change and coefficient of thermal expansion only. This assumption is 
based on the material property data for typical ablation materials 
which are reported in References 9 snd 10. These data indicate that 
the ablation materials have Young's moduli which are much smaller than 
those of conventional structural materials such as steel, or aluminum. 
Therefore, the substructure is much stiffer than the ablation material 
and stresses in the ablation material will not produce significant 
substructure displacements. The exterior surfaces of the ablation 
material are assumed free of external loads. This assumption is 
realistic because the aerodynamic forces on the heat-shield surface 
are negligible during space flight. 

The analysis assumes that temperature variations occur only 
through the thickness of the ablation material, the z-direction. This 
assumption is also reasonable because heat-shield cooling results from 
radiation to deep space from the exterior x-y surface of the ablation 
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material (see Fig. 3) and the cooling rate will be uniform over this 
surface. 

The ablation materials considered in this analysis are uniform 
density mixtures without honeycomb reinforcement. These materials are 
assumed to be isotropic and elastic with the properties E and a 
being functions of temperature only. Except for the initial cooling 
period, the ablation material temperature distribution will be fairly 
independent of time. Therefore, a steady-state temperature distribu- 
tion is assumed. Heat-transfer analyses for calculating the tempera- 
ture distribution are available (Ref. 11, for example). 

For this analysis, it is assumed that the deformations of the 

x-z plane, which is located well away from the heat-shield corners, are 
described by plane strain theory. The ablation material strain in the 

y-direction (e-y) is assumed to be equal to “the substructure strain. 
Since the substructure is assumed free to expand or contract, the value 
for e y is equal to (Xg&Pg. 

2.2 Governing Equations 

Figure 4 shows the mixed boundary conditions on the x-z plane 
which are analyzed for u and w displacements. The bonded edge of 
the ablation material -b < x < + b, z = O' has prescribed displace- 
ments while the other boundaries have prescribed stresses. Since the 
temperature and material properties are constant in the x-direction 





u = xa AT 
s s 

w - 0 


assumed, boundary conditions 
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and boundary conditions on the x-faces are identical, the z-axis is a 
line of symmetry. These symmetry conditions are 

u = u 2 = 0 on x = 0, 0<z<h (l) 

Therefore only the rectangle 0 < x '< b, 0 < z < h; is considered. . 

The analysis is formulated in terms of displacements rather than stress 
functions because displacement boundary conditions can be prescribed on 
the x-axis, but stresses on this axis can only be determined from the 
thermal- stress solution. 

From Reference 12, the two-dimensional stress equilibrium 
equations are 


A(* ) + A( T ) 

dx x Sz xz 


X(o ) + A(t ) 

5z z Sx zx 


= 0 


= 0 


( 2 ) 


The thermoelastic relations between stress and strain are 
(Ref. 12) 


c? x = A® + 2Ge^. 


ccE £1 
1 - 2v 


Cy = Ae + 2Ge y 


rr.F AT 
1 - 2v 


( 3 ) 


c z = Ae + 2Ge - 
z z 1 - 2v 

T xy = G Axy> T yz = G 7yz> T xz = G ?xz 



15 


where 


A = 


vE 


(1 + v)(l - 2v) 


> G = 


E 


2(1 + v) 


and e = £ x + 6 y + € z 


Assuming linear strain displacement relationships, that is. 


c _ Su _ 3w • _ /Su Swi 

x "S’ z " S’ 7xz " Is SJ 


and substituting equations (3) into equations (2), the two-dimensional 
equilibrium equations for plane strain in terms of displacements may be 
written as 


5 

(1 - V 

)Eu x 

M 


VEVy 

, vEw z 

3x 

1(1 + v ) ( 1 - 

■ 2v), (l + 

v ) (l - 2v) 

(l + v)(l - 2v) 


Ea AT 

+ 

d 

Eu z 


4- _ 

Ew x 

= 0 

(1 - 2v) 


Sz 

2(1 + 

v) 

2(1 + v) 

S' 

Eu z 

- + 

Ew x 

4. 

A 

(l - v)Ew z 

chc 

2(1 + v ) 


2(] 

L + v) 



(1 + v ) ( 1 - 2v) 


; vEUy vEvy Eg AT =Q 

(l + v)('l - 2v) (l + v)(l - 2v) 1 - 2v 

The term Vy in equation (4) is retained because the assumption ey = 0 
is more restrictive than necessary for plane strain (Ref. 13). . 

Therefore, provision for retraining ey as an assumed value of 



substructure strain only, is included in the equilibrium and stress 
equations for plane strain. 

Since the ablation material is rigidly bonded to the sub- 
structure, the boundary conditions at the interface are 

u - cc, s £!r s x, w = 0 on 0 < x < b, z = 0 ( 5 ) 


The stress boundary conditions yield 


(l - v)u x + vey + vw z - (l + v)a AP = o| x = b 


u z + w x = 0 


0 < z < h 


( 6 ) 


vu x + V€ 
U Z + w x 


y 

= 0 


+ (l - v)w z - (l + v)a ZbT = 0 


The boundary conditions for the corners at x = 0, 
z = h are 


0 < x < b 
z = h 

z = h, and 


(7) 


x = b. 


u = 0 

vu x + V£y + (1 - v)w z 

U z + Wjr = 0 




(l + v)a = 0 > 


x = 0 
z = h 


J 


( 8 ) 


(1 - v)u x + VSy + VW Z - 
VU^ + V£y +(l-v)w z - 
u z + w x = 0 


,(l + v)a £T = 0 
(1 + v)a = 6 


x = b 
z = h • 


j 

Equations (4) are integrated by finite differences.. 


( 9 ) 


Thus the 


x-z plane is subdivided by a rectilinear mesh (see Fig. 5) with a 
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constant mesh width Ax in the x-direction and a constant mesh width 
Az in the z-direction. Mesh spacing with Ax ^ Az is used because 
the structures being analyzed typically Tiave much larger dimensions in 
the x-direction than in the z-direction. The mesh lines are parallel 
to the coordinate axes with mesh lines on the edges x = b and z = h. 
The digital computer used, to solve the finite difference , 

equations did not permit the use ,of zero as a subscript. Therefore, 

♦ 4 > 

the mesh points x m and z n are defined by 

= (m “ 1 m — 2^ ♦ • ♦ f z^ — (n — X)^z ^ n — X^ 2^ • • • y J 

(10) 

■2 . 3 Finite Difference Form of Equilibrium Equations 
The method used to integrate- the differential equations of 
equilibrium by means of finite differences is given in Appendix A. 

A digital computer program for solving the finite difference equations 
was witten as part of this investigation. This computer program is 
discussed in Appendix B. A summary of the finite difference equations 
is presented in this section. Refer to Figure 5 for locations. For 
programing the finite difference equations for computer solution, the 
designations A, B, C, and D are used for recurring combinations 
of Poisson's ratio. Although this analysis assumes plane strain, plane 
stress solutions may be obtained from the computer program by using 
different values for A, B, C, and D. These constants for both 
plane strain and plane stress are given in Table I of Appendix A. 
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2.3-1 Corner Point 1,1 


u = 0 

(ID 

v = 0 

2.3-2 Axis of Symmetry m = 1, 2 < n < J - 1 

u = 0 


A(AEw z + BEe y + BEujj- - CEa AT) + (DEu z + DEw x )* = 


A 

'(En + E n+1 ) 

At ?■ 

2 


/ \ (^n ®n-l)/ 

( w m,n+l ” w m J n-l)-“ o (%,n “ w m : 




+ PA-^ Av. [?n + l u in+l^ n + l ” ^n-l^jn-i-l^n- J + 7^7 (^n+l " ®n-l 3 

" Rto zXE )n+l " ( s<x ^)n-ll + o?^T ~ f%i+l,n+l " %i+l.,n-ll 

PA v. ^ -J PAx; Az ^ 

BEy, »— — I 

+ 2 + ^rn+l, rx] = 0 

Ax^ 

2.3-3 Corner Point m = 1, n = J 


( 12 ) 


,u = 0 
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A(AEw z + BEe y + BEu^ - CEcc, ASC)' = 


A 


3Az c 


9(E n + E n-.i) ( ,_ _ ^ , (E n -1 + E n _2-)^_ 

2 (%,,n " ™m,n-l) + 2 ( w m,n-l ” w m, i 


+ 


B 


6Ax Az 


Q7Eri%.+l,n “ + 3En-2 u m+l_, n-2^j 


Bei 


+ g[-7®n - UZn-1 + 3Bn- 2 ] - ^['7(0* 


- 12 (Ea AT ) n _ 1 + 3(Ea Aff) n _ 2 ]] = 0 


2.3.4 Restrained Boundary 2 < m < I, n = l 


u = a s Ar s x 


w = 0 


2.3-5 Interior 2 < m < I - 1, 2 < .n < J - 1 


^-(AEu x + EE<= V + BEw z - CEa AT) + ^-(DEu z + DEw x ) 
ax J Acr 


8z 


AE n _ BE^ 

— 2 [%i--l, n ~ ^m^n + ^m+l^nj + 4 ^. \ y , [%i+l } n +1 “ w m+l,n-l 


“ w m-l,n+l + w m-l,n-l] + 4 ^ ^ 0 n+l^n+l,n+l ” w m-i,n+l^ 


®n-l( w m+l,n-l ~ w m-l_,n-lQ + ‘ ; 


®n+l) 

2 — v-n m ^ n +x 


- Vn> - ^(\n - Vn : 'i')] - 0 ' 


( 13 ) 


’ (l4a) 
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^-(AE-w z •+ EEu x + BE€y - CEa AT) + ^(DEu^; + DEw x ) = 
dz dx 


"(En+l+En),.. .. , (E n + En-l) 

. - u%,n+l “ “ 2 


.,n ’ %i,n-l)J 


4Ax Az 


(En+l(%i+l,n+l “ %-ljn+l) En-l (■%+!, n-1 “ %-l,n-l)3 


+ Sfe+HL " ^-l] " - (Ea AP) n _i] 


4Ax Az 


(%+l^n+l “ %-l,n+l “ %+!, n-1 + ’ u m-l,n-l3 


UJiiyi f** ’ 

+ p [%-l,n ~ ^V,n + ^%-l,nJ = ® 
Ax 


2.5*6 Exterior Boundary 2 < m < I - 1, n = J 


(lift) 


^-(AEu x + BEe T + BEw z - CEa AX) + ^-(DEu z + DEw x ) = 
dx J . dz > ’ 

( n B 2 \® n r ^ , i , D T 9(E n + E n . 

\ " A fe^ 1 " 1,n " ^2 1 • ' 2 


v ; 


fv, _ „ ) + (E n -1 + E n _ 

^n-ll 2 


^^“^,n-l " V>n-2) ' 

r , ( ' t , 


+ 1pAv ~G'^n( , %+i,n “ ™m-l,n) ” 9E^_x( w m+l, n -lf " w m-l, n-1 ) * • 

+ 3En_2 ( w m+l, n -2 “ %-l,n-2)] = 0 , ' (^a) 
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J^(AEw z + BEu^. + BEe y - CEa Zffl) = 


A 


JAz £ 


9(E n + E n _]_) ^ ^ (En-1 +E n-2),\ r , T - T . T 

2 K,n " ^V,n-1' 2 '^Vjii-i “ ^ J n-2'' 


B 

r 

12Ax 



^•TEn (%+X,n “ %L-l^n) “ (%+l, n-1 “ 

i 

+ 5E n _2(% 1 +i ,n-2 “ ‘ u m-l,n-2)3 + ^— •’(^7En “ ^®n-l + ^n- 2 ] 

- ^-[^7(Ea AE) n - 12 (Ea ^P) n _i + 3(Ea £T) n _ 2 ] . = 0 , * ' - (l^b) 


2.3.7 Exterior Boundary m = I, 2 <j n ^ J - 1 


g— (AE u^. + BESy + BEw z - CEa £T) - 


+ - 1 - ( ^ u m-l ? n “ V1 m-2 > n3 + 1 pA-y a ?. E ^ n+1 “ 


3Ax‘ 


~ 12 ( w m-l,ri+l “ w m-l,n-l) + 5 '( w m-2,n+l ~ v m-2,n-l0: 


8 BEn€ y + 3 ^ c(Ea ^ = 0 


3Ax 


(l6a) 



B 2\ 1 

(®n + ®n+l) 

A L 2 

2 


/ \ (^n + ®n-l) r 

v"%,n+l " 2 " ^^n-l/ 


DE DE 

+ + ■ LD %-l,n ~ %-2^n3 + Vpa^ a?' (^ n+l " ' u m.,n-l) 

12(%-l, n +l - %-l,n-l) + 3(%~2,n+l “ "%-2,n-l)] 


‘ ‘ f )&*" ffl W - (®» ffl >n-l] 


+ £ 1 - !')&*«• - •*£ - ° 


(l6t>) 


2.3*8 Corner Point m = I, n = J 


— (AEu x + EEey + BEw z - CEa AT) = 


( A ’ t'^n + 1 ° Um - 1 > n " ^2,n] 

+ _8_f A . *\r g«fy * C(Ect „ o 

3Ax\ A /{J A + B J 


(17a) 
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^-(AEw z + BEu x + BEey - CEct £3?) = 


A - 


B‘ 


(En-1 + 


>\ x [ 9 g^JM^) ( . „ * ) + 

/ 3Az 2 L 2 . m ' n m,n ~ 1 2 


( w m,n-l " w m,n-2^ 1 



jg^ “ 4E n-l + v 2 


+ C (f ‘ «>n -'M&i ^>n-l. + 


♦ -§-/* 

3Az\ A 


a/L 


BE- 


l n e y + C(Ea £E) n 


A + B 


= '0 


®n-2 ) 


(17b) 



CHAPTER III 


EXPERIMENTAL TEST PROGRAM 

A test program of limited scope was conducted to investigate the 
behavior of ablation materials subjected to restrained shrinkage. 

Eight small models consisting of rectangular panels of ablation 
material bonded to steel substructures were tested in a cryochamber 
at the Langley Research Center. 

3*1 Model Configuration 

The model configuration is shown in Figure 6. The models were 
fabricated by bonding two k- by 11-inch ablation material panels to a 
O.OSO-inch-thiek steel plate. The steel plate was made wider than the 
ablation material panels to provide a means to mount the model in the 
cryochamber. The planform area of the model was necessarily much 
smaller than an actual heat shield, but the ablation material thickness 
was full scale. A symmetrical model was used to eliminate bending. 

The model substructure was sandwiched between two ablation material 
panels which were both cooled simultaneously. The model substructure 
provided an extensional stiffness similar to an actual structure. 

This model configuration furnished a shrinkage restraint boundary 
condition on one surface of the ablation material and stress free 
boundaries on all other surfaces. These boundary conditions 
simulate those of an actual heat shield and agree with those 
assumed in the analysis. To indicate ablation material cracking 
during the tests, crack detection circuits of conducting paint 
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Figure 6.— Test model configuration 
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were applied to the model surfaces. These crack detectors are visible 
in Figure 6. The models were also instrumented with thermocouples to 
measure the temperature of the ablation material and substructure 
during testing. Neither the thin layer of conducting paint nor the 
thermocouples bonded into the ablation material appeared to have any 
physical influence on the cracking of the ablation material during the 
tests. 


3.2 Ablation Material 

Six models were fabricated with the ablation material shown in 
Figure 6. This material is a low density composite of nylon and 
phenolic which has been widely tested as an ablation material (Ref. 14) 
and is typical of the materials considered for ablative panels on 
future manned spacecraft. The phenolic nylon material was tested at 
thicknesses of l/% 5 / 8 , and 1 inch. Three models were bonded with an 
adhesive having a cure' temperature of 325° F and the other three were 
bonded with a room temperature curing adhesive. Both adhesives are 
commonly used in heat- shield fabrication, particularly the high 
temperature curing adhesive. The thermophysical properties of the 
phenolic nylon test material are reported in Reference 10. Data from 

Reference 10 ~ which are pertinent to this investigation are shown in 

■» 1 * ' 

► < 

Figures 7 a* 1 *! 8. The faired curves shown in these figures are also 
from Reference 10. These data show large variations and indicate 
uncertainties in derived strength properties. Such variations in data 

-J-These data are presented with the written permission of the 
authors. 
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Figure 7*— Elastic modulus and Poisson’s ratio versus temperature for 
phenolic— nylon ablation material loaded in tension 
(data from ref. 10 ) 


Ultimate tensile strength, psi 


27 



Temperature — degrees F 

Figure 8.— Ultimate tensile strength versus temperature for phenolic- 
nylon ablation material (data from ref. 10) 
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are typical for ablation materials. Phenolic nylon ablation material 
also exhibits dimensional instability when subjected to changes in 
humidity. The constituent materials are hydroscopic and will absorb 
atmospheric moisture which produces swelling. Exposure to vacuum will 
cause the material to lose moisture and shrink. 

Of the remaining two models, one was fabricated with NASA 602 
and one with the Apollo spacecraft ablation material. NASA 602 is a 
silicone elastomer material described in Reference 1. This material 
differs from phenolic nylon in composition and is less brittle at room 
temperature, but has a coefficient of thermal expansion which is about 
twice as large as phenolic nylon. A 1-inch thickness of NASA 602 was 
tested and it was bonded to the substructure by using the constituent 
silicone resin as an adhesive. The Apollo heat-shield material is an 
epoxy -based composite reinforced with a plastic honeycomb. Heat 
shields of this material are fabricated by bonding the honeycomb to the 
spacecraft structure and then filling the honeycomb cells with ablation 
material. The resulting ablative shield consists of many small columns 
of ablation material separated by honeycomb cell walls. A 1-inch 
thickness of the Apollo material was tested to compare its behavior to 
the phenolic nylon. 


3*5 Test Procedures 

The models were tested in a cylindrical vacuum chamber whose 
walls contained liquid nitrogen. The models were supported in the 
test chamber so as to permit simultaneous radiation cooling of. both 
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rectangular slats of ablation material. The test procedure was as 
follows: The model was placed in the test chamber and the instrumenta- 

tion was connected to recorders. Radiation shields consisting of 
fibrous insulation with an aluminum foil cover were placed over all 
exposed surfaces of the steel substructure and the model ends. Only 
the exterior surface of the two ablation material slabs was exposed 
to the cold chamber walls. The test chamber was sealed^ evacuated to 
a pressure of about 1 torr, and liquid nitrogen was then introduced 
into the chamber walls. 



CHAPTER IV 


RESULTS AND DISCUSSION 

The test program described in Chapter III was performed to 
obtain data from typical ablation materials which would identify the 
critical stress regions during restrained shrinkage and which could be 
qualitatively compared with the analytical stress analysis. The 
computer program developed for this investigation was used to calculate 
stresses caused by restrained shrinkage. These calculations have been 
compared with exact solutions. Comparisons between the calculations 
and the test data have also been made. The results of these studies 
are presented and discussed in the following sections. 

4.1 Experimental Results 

The phenolic nylon plate models described in Chapter III are 
shown after testing in Figure 9- The models shown in Figures 9(a) and 
9(b) had the steel substructure sandwiched between two 1/4- inch thick- 
nesses of ablation material. These models developed cracks at a surface 
temperature of -125° F and substructure temperature of -120° F. Thermal 
stress calculations, which are presented in a subsequent section, 
indicate that this cracking occurred in the region of maximum in-plane 

tension stress. Edge cracking was not observed. 

Figure 9(c) shows the interface cracks across the narrow 

end of the 5/8-inch thick models. These models cracked at a 
surface temperature of 30° F and substructure temperature of 72° F. 
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(a) Ablation material thickness of l/4 in., bonded with 
room temperature curing adhesive 

Figure 9 .— Post test photographs of phenolic— nylon plate models 
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(b) Ablatioii material thickness of l/4 in., bonded with 
high temperature curing adhesive 


Figure 9*“ Continued 





(c) Ablation material thickness of 5/8 in., bonded with 
‘ room temperature curing adhesive 


Figure 9 *~ Continued 
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(d) Ablation material thickness of 1 in., bonded "with 
high temperature adhesive 

Figure 9 .— Concluded 
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The crack location indicates that the cracking was caused by edge 
tension stresses. The cracked region of the models having 1- inch- thick 
ablation material is shown in Figure y(d) . The crack circuits indi- 
cated that these models developed end cracks while the test chamber 
was being evacuated prior to cooling. As discussed in Chapter III, 
phenolic nylon shrinks from moisture loss during vacuum exposure and 
this shrinkage produced edge stresses large enough to crack the 
ablation material. Subsequent cooling of the model to a surface 
temperature of -75° F increased the separation and caused the large 
crack visible in Figure 9(d). 

The models with the ablation material bonded with high tempera- 
ture curing adhesive performed essentially the ’Same as those fabricated 
with a room temperature curing adhesive. The bonding operation at 
32 5° F apparently did not cause significant thermal stress in the 
models. The reason for this is that the strength of this phenolic 
nylon material is extremely small at 525° F (see Figs. 7 and 8). 
Therefore, it is anticipated that the hot ablation material accommodated 
the differential contraction strains and did not develop significant 
thermal stresses during subsequent cooling to room temperature. 

The test model fabricated with 1-inch- thick NASA 602 also 
developed end cracks at the location shown in 'Figure 9(d). The crack 
circuits did not function well with this material and, therefore, the 
temperature at cracking was not determined. This test result was 
expected and indicates that silicone elastomers are also subject to 
cracking from restrained shrinkage at low temperature. The model 
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fabricated with 1- inch-thick Apollo heat-shield material did not 
develop cracks during testing. This result was also expected because 
the honeycomb reinforcement is effective in preventing cracking due to 
restrained shrinkage. 


4.2 Analysis Results 

The accuracy of the numerical analysis is evaluated by comparing 
the calculated stresses with an exact solution. The exact case is the 
plane strain analysis of a semi-infinite slab with constant material 
properties. The slab is subjected to a uniform temperature change and 
is restrained in both the x- and y-directions such that the strain in 
these directions equals a g AT s . From Reference 15, the solution for 
this case (away from the edges) is a x = by = - Ea. AT/l - v, >o z = 0 
where a is the difference between the slab and substructure expansion 
coefficients. — - 

Calculated stresses in the slab subjected to a uniform temper- 
ature decrease are shown in Figure 10. The stresses are normalized with 
respect to the parameter Ea £T/l - v and the abscissa is the ratio of 
distance from the free edge to the slab thickness. As the distance from 
the free edge increases, the edge effects diminish and the calculated 
stresses agree with the exact solution in regions for which the exact 
solution is valid. Mesh spacing s of Ax = 0.125 in. and Az = 0.025 in. 
were used for all the calculations presented in this thesis. 

The slab thermal stresses at the interface z = 0 are shown in 
Figure 10(a) and serious discrepancies exist between the calculated 
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(a) Stress distribution at ablator— substructure interface, z - 0 


Figure 10.— Calculated restrained shrinkage stresses for a uniform 
temperature slab 
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Stress distribution at the slab mid-depth, z 
Figure 10.— Continued 


39 



(c) Stress distribution at the slab surface, z = h 
Figure 10.— Concluded 


values and the free edge stress "boundary conditions,, o x = t xz - 0, at 
x - "b. This result was not anticipated and the analysis was completed 
before the results from the displacement analysis of Reference 8 were 
known. Much effort was spent in trying to resolve the difficulty, 
including the use of extremely small mesh spacing s, but no successful 
solution was obtained. It is noted that corners are excluded in 
studies of finite difference approximations for mixed boundary-value 
problems (Refs. 16 and 17, for example). A stress singularity at the 
corner point is indicated by the large shear stress values obtained 
with the present analysis. The analysis of Reference 18 also identifies 
the corner point as a potential singularity point. 

A stress singularity at the comer point, x = b, z = 0, means 
that the stresses at that point cannot be evaluated, but it does not 
invalidate the results obtained at other regions of the slab. It is 
emphasized that the singularity is one of stress, but not of displace- 
ments. With the present analysis,, the interface displacements are 
specified as u = a s ,AT s x and w = 0. These values for u and w 
also specify the values for u^. and w x . The interface displacements 
appear as known values in the equilibrium equations for the mesh points 
at 1 < m < I, n = 2 and these displacements, therefore, are used in 
the systems of equations which are solved for u and w throughout 
the x-z plane. Since the interface displacements are known, however, 
the equilibrium equations are not applied at the interface. This 
means that the prescribed stresses at the comer point are not included 
in the slab displacement equations. The interface stresses are 
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calcxilated -with, the specified values for u^. and w x and from values 
for u z and w z obtained with forward difference equations using the 
specified interface displacements and the calculated slab displace- 
ments. Since w = w x = 0 is specified along the interface, satisfying 
the shear stress boundary condition requires that u z = 0 at x = b, 
but this condition is not in the slab displacement equations. Since 
the slab contracts more than the restraining substructure, the u 
displacements at the mesh points 1,2 and 1,3 are greater than 
Uj The u displacements at these three points do not lie along a 
parabolic curve. Therefore, a three point forward difference equation 
gives a negative and nonzero u z . Decreasing the mesh spacing, Az, 
merely increases the calculated value for u z at x = b. To satisfy 
the cr x = 0 condition, with u^. specified, requires negative values 
for w z which are considerably greater than the values at x less 
than b. To satisfy the slab equilibrium equations along the x = b 
edge of the x-z plane, however, w x must be positive at the mesh 
points for n > 2 because u z is negative at these points. Therefore, 
the forward difference equation does not give large negative values 
for w z at the corner point. The net result is that g z at the 
corner point is a tension stress rather than the compression stress 
given by the free edge boundary conditions. 

It is noted that for slab equilibrium, the sum of cr z stresses 
along the interface must be zero. The o z stresses shown in 
Figure 10(a) satisfy this requirement, however, this does not indicate 
a unique solution. It is speculated, however, that even if compressive 
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a z stresses exist at the comer point, a stress reversal must occur 
near the corner. In Reference 5 the author states that photo elasticity- 
measurements up to 0.001 inch from the comer showed tension values for 
cr . The test results from the present investigation also indicate 
tension stresses as evidenced by the edge cracking at the interface. 

The numerical plane- strain analysis gives stress distributions 
in the x-z plane which: (l) agree with the exact solution well away 

from the edges; (2) agree with the prescribed boundary stresses at the 
slab mid-depth and surface (Figs. 10(b) and 10(c)); and (3) qualita- 
tively agree with the plane stress analyses of References 4 and 8 at 
the mid-depth and surface locations. In Reference 4, the author states 
that the effects of the stress singularity were highly local and that 
valid solutions were obtained away from the singularity. Because of 
(3) this statement can be applied to the present analysis. Also, the 
displacement equations are solved by repeated iterations of simul- 
taneous equations and if significant perturbations were introduced by 
the comer effect then these perturbations would influence all displace- 
ments. Such perturbations would cause the calculated stresses in all 
regions of the x-z plane to diverge from exact solutions. Since this 
divergence does not occur, it is concluded that the comer effects are 
confined to the immediate vicinity of the comer and that valid stresses 
are calculated at other regions. 

Figure 10(b) shows that at the mid- depth, z = h/2, of the 
restrained slab the edge stresses are minimal. Therefore, to investi- 
gate edge tension stresses which might cause the ablation material to 



crack., it was necessary to select an edge region which, was somewhat 
removed from the comer point, x = b, z = 0, and yet in the highly 
stressed region. A location at z = 0.05 inch was selected. At this 
location, good agreement was obtained with the prescribed boundary 
condition that a x - 0. The calculated shear stresses at this boundary 
location reached peak values near the edge and indicated a rapid 
decrease toward the prescribed edge value. Figure 11 shows the 
calculated values for c z , normalized with respect to Ea AT/l - v, 
for various thickness slabs subjected to a uniform temperature 
decrease. These values are calculated at x = b, z = 0.05 inch, and 
show that the edge tension stress increases with slab thickness. 

This result agrees with the experimental data presented herein. 

In actual applications, an ablation heat shield would have a 
temperature gradient in the z-direction rather than the uniform 
temperature assumed in Figure 11. Therefore, the influence of 
three different temperature distributions on restrained thermal 
stresses was investigated. Case I assumed a 'uniform temperature 
distribution. Case II assumed a linear temperature distribution through 
the ablation material thickness with AP = 0 at z = 0 and Case III 
assumed a uniform temperature in the upper half -of the slab thickness 
(h/2 < z < h) and a parabolically increasing temperature to AT = 0 
at z = 0 in the lower half of the slab (0 < z < h/2). The value for 
AP at z = h was the same for all three cases. The temperature 
distributions for Case II and Case III are suggested by the temperature 
distribution measured during model tests. The calculated stresses are 



(b) h = O.35 in. 


Figure 11.— Calculated a stresses for various thickness slabs 
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at z = 0.05 in. (slab subjected to restrained shrinkage 
at a uniform temperature) 
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compared in Figure 12 -which shows the ratio of a z at x = b, 
z =0.05 inch to crx max f° r various assumed slab thicknesses. 

Figure 12 shows that the edge tension stress increases with slab 
thickness for all three temperature distributions. The Case III 
temperature distribution induces the greatest slab curvature and 
causes the highest edge tension stress near the nonbending interface. 

As expected, the linear temperature distribution, of Case II produced 
the lowest tension stress. The' most interesting point. about these 
calculations is the indication that the region of maximum .tension 

f 

stress shifts from the slab interior to the edge when the slab thick- 
ness exceeds about 0.5 inch. • The calculations are consistent with the 
test results. Interior cracks developed in the test models with a 
slab thickness of l/4 inch and Figure 12 shows a value of o z ) a x 

ni8iX 

less than 1.0. The models with slab thickness of 5/8 and 1 inch 

developed edge cracks and the value of for these thicknesses 

z x max 

is greater than 1.0. 

4.3 Comparison of Calculations and Test Results 
Thermal stresses were calculated for a phenolic nylon material 
at the temperature when model cracking occurred. The test data 
indicated that the model with l/4-inch- thick ablation material 
(Fig. 9(a)) cracked when the surface temperature was -125° F and the 
substructure temperature was -120° F. The reference temperature was 
80°. The material properties were obtained from data in Reference 10. 
These data indicate that at the measured temperature the thermal 
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expansion coefficient is 23-5 X 10 “6 in/in/deg F. However, as shown in 
Figure 7 the data for Poisson’s ratio and Young ' s * modulus show wide 
variations. The values selected were: v = 0.1 and E = 120,000 psi. 

The expansion coefficient for the steel substructure was taken as 

S 4 

6.76 x 10“° in/in/deg F. The thermal stress calculations indicated a 
maximum tension stress of 460 psi at the top surface of the ablation 
material. A calculation assuming Poisson's ratio equal to 0.25 
produced a maximum stress of 550 psi. Examination of the ultimate 
tensile strength data in Figure 8 shows that these stress values would 
cause ablation material cracking such as occurred in the model tests. 

Thermal stress calculations were also made for comparison with 
the 5 / 8 -inch model (Fig. 9(c)). The model cracked at a surface 
temperature of 30° F, a mid-depth temperature of 32° F, and a sub- 
structure temperature of 72° IT. The reference temperature was 80° F. 
The material properties were the same as previously discussed with the 
exception of Poisson’s ratio which was assumed to be O. 15 . The maximum 
tension stress occurred at the edge, but was only I 50 psi, and well 
below the ultimate tension value from Figure 8 . Apparently the stress 
resulting from moisture-loss shrinkage was primarily responsible for 1 
the model cracking and not the stress from the relatively small 
decrease in temperature. 



CHAPTER V 


CONCLUDING REMARKS 

An approximate numerical analysis for 'restrained, shrinkage ‘ 

stresses in rectangular slabs has been presented. Test results from 

J % 

rectangular slabs of ablation material subjected to severe restra in ed 
shrinkage have also been presented and qualitatively compared with 
analysis calculations. The present analysis which is programed for 
solution with a digital computer calculated restrained shrinkage 
stresses which are consistent with the test results reported herein 
and with other experimental investigations. The alternating direction 
technique (Ref. 19) which was programed to solve the finite difference 
equations appears to be rapid and efficient. 

The analysis was an attempt to solve the mixed boundary- value 
problem of a rectangular plane having both specified displacement and 
stress boundary conditions. The analysis objective was not fully 
achieved because no satisfactory solutions were obtained at the 
boundary corner point for which both displacements and stresses were 
specified. The numerical analysis did not determine the strains and 
stresses at this boundary point with any precision. Extensive and 
time consuming efforts were made to improve the corner point solution, 
but these efforts were unsuccessful. It is not known what modifica- 
tions to the analysis would be required to obtain better solutions. 

Both the test results and the stress calculations from this 
investigation indicate that a heat shield subjected to restrained 
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shrinkage can develop edge separations between the ablation material 
and the substructure. Heat shields fabricated with bonded, 
unreinforced ablation materials appear susceptible to edge separation 
or cracking when the ablation material thickness exceeds about 
0.5 inch. The data also indicate that this cracking may occur merely 
because of changes in the moisture content of the ablation material 
and at temperatures much less severe than those which could be 
encountered in space flight. The desirability of using honeycomb 
reinforced ablation material to withstand restrained shrinkage is 
evident, but must be weighed against fabrication cost and difficulty. 

A reinforced region at the heat-shield periphery might be an acceptable 
compromise solution. 
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APPENDIX A 


FINITE DIFFERENCE EQUATIONS 

In this section, the partial differential equations of equili- 
brium are expressed in finite difference format suitable for numerical 
calculations. The resulting series of difference equations is solved 
by an iterative procedure which is discussed in Appendix B. 

The finite difference expressions for second derivatives are 
obtained using the "half station" approximation which is the natural 
result of making the finite difference approximation before expanding 
the derivatives. This approximation is used because the modulus of 
elasticity is variable and the results of Reference 20 indicate that 
the "half station" approximation is more accurate than the "full 
station" approximation. Accordingly, 


dz (EWz) m,n Az ^ z ^m,n+ l/2 " m,n- 1/2 

= ^2 j^ cl+ l/2 w m, n+1 " (®ri- l/2 + E n+ l/2^ w m,n 


+ E 


'n- l/2 w m,n-l 


(A-l) 


The values of E are evaluated at the temperature of a point midway 
between grid points m,n and m,n+l; m,n and m,n-l. Thus 
equation (A-l) becomes 
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The finite difference expressions for first derivatives are obtained 
with the "full station" approximation. Three-point central difference 
equations are used at interior mesh points and three- point forward 
differences are used at exterior boundary mesh points. 


Interior Points 

The equilibrium equations for plane strain are. 
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• (A- 3b) 


The values for A,B,C and D are given in Table I. 

The finite difference foim of these equations is obtained as follows: 
The first term in equation (A-3a) is obtained by a three-point central 
difference approximation 
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The second term of equation (A-A), a mixed derivative, is obtained by- 
making a full station approximation in the x-direction before expanding 
the derivative. 
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The third and fourth terms in equation (A-3a) are constant in the 
x-direction and, therefore, are dropped. The fifth tern in 
equation (A-A) is obtained as shown in equations (A-l) and (A-2). 
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v%i,n+l - u m,n) ( u m,n - %i,n-l) 


. (A-6) 

The sixth term in equation (A4) is obtained in the -same manner as 
equation (A-6) except that E is considered to vary in the z— direction. 

f 

— (DEv x )^ n = ^ — (jn +1 (w m+1#n+1 - w m _ ljn+ i) - E n-i^ w m +l,n-l " w m-l,n-lD 

(A-7) 


By combining equations (A- A) through (A-7) the finite difference form 
of equation (A-3a) becomes 



5T 


AE n . 


BE 

+ u m+l, n] + "4^ ^ {™m+l,n+l ~ w m-l,n+l 


/sc £ 


- w 


, -1 , D | ~( E n + E n+l) ( v 

,n~l m-l,n-lj ^ 2 \_ 2 Hi,n+1 ~ 


(Eji + E n _i) 


(«m 


t,n “ + ^2 (|n+l( w m+l,n+l “ 1 %-l,n+l) 


S n-l^ w m+l,n-l “ w m-l,n-lO ~ 0 ( A ~8) 

r similar procedures equation (A- 3b) becomes 

■(E n +E n+ i),__ „ , ( E n +E n-l)/ 

2 ' w m.n+l %,n^ 2 ^ w m,n " %, n-1 

+ 4 Ax'~Az 6 n+1 ^ +1 ^ n+ l “ %-ljfc+l) ~ E n-l( u m+l,n-l ~ ' u m-l > n-l)] 

+ - ®n-J - - <&■ «)n-l] 


DE. 


1 p ~i 

+ ]^L^l,n + l “ ‘ u m-l,n+l “ ^m+l^n-l + tl m-l J n-lJ 
+ ^gl^n-l^n ~ n + w m+l^ nj = 0 


(A-9) 


External Boundary 
2 < m < I - n = J 
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The following procedure is used to obtain the finite difference 
form of the equilibrium equations at the external boundary. The first 
term of equation (A-Ja) is 

-s 

— (AEu x “ ^n + %+l,n} (A-10) 

The second term in equation (A- 3a), b/bx (EEw„), is obtained using the 

Yt 

normal stress boundary condition, that is, 

cr = 0 = AEw + EEu + EEe - CEa (A-ll) 

z z x y ' ' 

This gives 


— (BEw z ) m ^ n = “ " g(%i-l,-n " ^%i,n + u m+l,n] (A- 12) 

The third and fourth terms in (A-3a) are constant in the x-direction 
and are dropped. The pattern for the fifth term b/bz (DEu z ) is 
obtained with a three-point "half station" forward difference 
approximation. 


^~(DEu z ) 
oz m, n 


= 8Eu, 
3Az 


2 - 9Eu z , + Eu 

z m, n z rn, n- 1/ 2 nn, n- 3/2 


(A-13) 


Applying the shear-stress boundary condition, (Eu z )^ ^ = -(Ew x )^ ^ 
equation (A-13) becomes 
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D 

3Az 


+ 9<En * g *»*W - %,n-l> 


+ (En -l * En ~ g) (%, n -1 - %,n-2) 


(A-14) 


The sixth term in (A- 3a) is obtained with a three-point "full station" 
forward difference approximation. 


a z (DEw x)m,n “ 4zh?Azl? En ^ %+1 ^ n " 

f 

^®n-l^%+l,n-l “ " H m-l,n-l) + E n-2^ w m+i J n-2 “ w m-l, n-2^ 

(A-15) 

By combining equations (A- 10% (A- 12% (A-l4% and (A-15) the finite 
difference form of the equilibrium equation (A- 3a) becomes 


5 2\E. 


A " A ]J2 n " ^ n n ] "77 


3Az £ 


9(E n +E n _!) 

2 (%i,n “ ' u m,n-l) 


’(El ) 

+ n - 1 ”(■%, n-1 “ %!, n-2 ) + tpav"" A 7 . E^ 51 ( Wm+ 1> n ~ w m-% n ) 


" ^-l(%+l,n- 1 _ %-l,n-l)‘ + 3Ba-2(%+l,n-2"' H ia-l,n-2j] = 0 


(A- 16) 


The finite difference form of the second equilibrium equation (A- 3b) 
is obtained in a similar- manner. This equation is 



6o 


A 9(E n + E n-1 ) ( E n-1 + ^- 2 )/ \~| 

2 ~ 2 '■ w in,n “ %,n-l) + “ '• w m,n-l " %,n-2J 


+ 1 pAv ^rt'^ n ^ tlni+ ^ n ” “ 22E n~l (^m+1, n-1 ~ ^ra-l^n-l) 

+ 5E n-2^%-fl,n-2 “ %-l, n-2 fl + " 12E n-l + 5E n-2] 

- ^-[-T(Ea /ff) n - 12(Ea £P ) n-1 + 3(Ea £T) n _ 2 ]' = 0 (A- 17 ) 


External Boundary 
m = I, 2 < n < J - 1 

The finite difference form of the equilibrium equations (A- 3a) 
and (A-5b) - OT this boundary is obtained in the same manner as 
equations (A-l6) and (A-17). The equations are 

+ 10u m-l,n " %-2,nj + tpa-v a* “ %n-l) 

" 12 (*5n-l,n+l “ 1 ia-l,n-l) + 5 ^-2,n+l “ V-2,n~l0 

' 3^ B ^* e J r + jlr c(Ea = 0 (All8) 
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B 2 \ 1 

( E n + ^+1). 

A/Az2 

2 


( w m,n+l “ w m,n^ 


(E. 


n + ^^( Wl 


m,n 


- w 


Wjn-1^ 


+ En 2 [;9w m ^ n + 10w m-1#n - w m _ 2 ,iQ + 07 ( n m, n+1 “ u m,n-l) 

12(u m _i^ n +i “ ^-l^n-l) + 3(iiia-2,n+l “ u m-2,n~l3 “ ?/\J (^ ~ ~p^j 
OEa ffl) n+1 - (Ex ar) n _ J + |g(l - ^ )g n+1 - Vl} * 0 ' (A-19.) ' 


Comer Point 


m = I, n = J 


At this comer point, the normal stresses a x , a z , and the 
shear stresses are zero. The equilibrium equations are converted to 
finite difference forms with the following procedure: With E and T 

constant in the x-direction and with zero shear stress the equilibrium 
equation (A-3a) becomes 


^•(AEu x + BEw z ) = 0 (A-20 ) 

dx 

After substituting the value for BEw z obtained from the boundary 
condition a z = 0, equation (A-20) is 



(A- 21) 
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Equation (A-21) is expanded with a three-point "half station" 
difference approximation. 


_ 6 _ 

5x 


k A -fK 


== A - 


m.n 


A/ 3 Ax 


8(u x ) -9(Ux) + (uJ , 

A m,n x m~ 1/2, n x m- 3/2, n 


(A-22) 


From the normal stress conditions at the comer point. 


( u x) 


EE n e y + C(Ea AT) n 


m,n 


A + B 


(A-23) 


The finite difference form for equation (A-21) is obtained by com- 
bining equations (A-22) and (A-23) 



+ 




BE n €y + 'C(Ecc, AT) n 
A + B 


(A-24) 


Similarly, the finite difference form of equation (A-3b) is 
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L - *U_ 

l a /3Az 2 


- 9(En * En ~ l) (Vn - ’Sn-l) + 


- %,n-2> - ( B - x)l&X ' + Vj • 

■ + °(f ■ ^>n - 4 < E “ “>n-l + ")n- 2 ] 


8 L 

- —'I 

BE n €y + C(Ea AT) n 

3Az\ 

A/ 

A + B 


= 0 


(A-25) 


Internal Boundary 

At the internal boundary the u displacement is zero. 
Therefore, only the w displacements are unknown and these' are 
obtained from the finite difference form of equation (A- 3b) which 
considering symmetry is 


A 

(E n + E n+1 ) 

Az 2 
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/- \ (®n + ®n-l)/ \ 

n+1 “ w m,n-l) g " ( v m,n ~ w m,n-l) 
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2Ax Az 


(j'n+l u m+l, n+1 “ E n-1%+1, n-J + T^fei+l " E n-l] 


C 

2Az 


DB r 


0 Ea AI) n+1 - (Ea + fem+1, n+1 " %+l,n-l] 


+ ”"§G 2 %,n + 2%+l,n] = 0 


(A-2 6) 
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Internal Corner 

The finite difference form of equation (A-3b) for this corner 
point is obtained in the same manner as equation (A-17) and is 
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3Az £ 


9(% + En-l)/ ^ , (En-i + En.g), ^ 

— < w m,n ” %n-l) + g <%,n-l " %n-2) 


B 


6Ax Az 


(j‘7 E n%+l,n " 12E n-l u m+l,n-l + ^ E n - 2 %+l,n -23 + 6az 


- 12^,! + 3^3] - g^[- 7 (Ea AT) n - 0^(Ea £SS) n -l 


+ 3(Ea AT) n _g] = 0 


(A-27) 


The finite difference equations given in this Appendix are 
solved for u and v displacements with a digital computer program 
which was written as part of this investigation. Some details of this 
computer program are' given in Appendix B. 























APPENDIX B' 


Tliis section describes the procedure used to solve the finite 
difference equations of Appendix A for the u and w displacements of 
each mesh point. The method for calculating stresses from these 
displacements is also discussed. 

Since the present analysis is concerned ■with edge stresses, it 
is necessary to use relatively small mesh spacing which in turn requires 
the solution of large sets of simultaneous equations. These equations 
are solved by an alternating direction method (see Ref. 19) programed 
for a Control Data Corp. 6600 Digital 'Computer. The alternating 
direction procedure is based on a line-by-line solution of relatively 
small sets of simultaneous equations that can be solved by direct, 
noniterative methods. With this method the finite difference formula- 
tions for equation (A-3a) are solved for- u displacements and 
similarly the w displacements are obtained from the finite difference 
formulations of equation (A-3b). >\ 

The displacement calculations are started by assuming values for 
the u and w displacements at each mesh point. The u and w 
displacements at mesh points in each horizontal line are calculated and 
stored. The displacements in each vertical line are then calculated 
and compared to the stored values. The .calculations of displacements 
are continued line-by-line until the absolute values from successive 
iterations agree within a specified percentage. The coefficient 
matrices for n and w m n are tridiagonal in both the horizontal 
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and vertical directions except on the boundaries x = b and z = h. 

At these boundaries the coefficients are modified by an algebraic 
transformation so that each set of simultaneous equations can be 
solved by a direct factoring method. 

The number of iterations required to achieve a selected con- 
vergence varies with the number of mesh points. A convergence 
accelerating factor is not used, but a procedure is programed whereby 
the displacements calculated for a coarse mesh spacing are' used as 
initial inputs for successively finer mesh spacings. This procedure 
produces significant savings in computer time, ho problems of 
stability are encountered since the program solves steady- state 
elliptic equations. 

'The normal and shear stresses at each mesh point are calculated 
with the finite difference forms of the following equations: 

, ♦ 

a x = + BEw z + BEe y - CEa AT (B-l) 

a z = AEw z + BEu^- + BEe y - CEa AT- (B-2) 

* r 

< i 

* 

t x Z , = D(,Eu z* + Ew J ' (®- 3 ) 

Equations (B-l), (B-2), and (B-3) are programed using three-point 
central differences for interior mesh points and three-point forward 
differences for external houndary mesh points. The computer program 
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uses the displacements computed for each mesh point to evaluate the 
strain values in equations (B-l), (B-2), and (B-3) with the exception 
of €y which is constant for this analysis. 



